Dynamics of correlations in shallow optical lattices 
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We explore how correlations evolve in a gas of lattice bosons when the lattice depth is rapidly 
reduced. We find a simple closed form expression for the static structure factor in the limit of 
vanishing interactions. The corresponding real-space density correlation function shows multiple 
spatial oscillations which linearly disperse as a function of time. By perturbatively including the 
effects of the interactions we calculate how the boson quasi-momentum evolves following the quench. 



Introduction. — While many phenomena in nature can 
be qualitatively understood by simple "mean field" -type 
theories, these approaches inevitably capture only a sub- 
set of the interesting physics. For example, the Mott in- 
sulating state of lattice bosons is not inert, rather there 
is a gas of particle-hole pairs which give rise to a finite 
correlation length. Similarly, a superfluid is not fully 
characterized by a coherent state. Recent experimental 
advances in ultra-cold gases have given us new tools for 
studying these fluctuations jH-Q. A particularly promis- 
ing technique (largely unique to cold atoms) is to rapidly 
change the Hamiltonian parameters such as hopping rate 
J and interaction strength U. The evolution following 
such a quench gives many insights into the single and 
many-particle properties of the system — the spectrum 
of excitations [9|, the manner in which correlations de- 
velop (lo| . and the role of quantum coherence Here 
we calculate how various correlation functions evolve af- 
ter a sudden quench from a strongly interacting Mott 
insulator (U 3> J) to a weakly interacting superfluid 
(J 3> U). Our calculations are closely related to ex- 
periments at Munich [l(|, but we consider a quench to 
much weaker interactions. 

We calculate how density-correlations evolve following 
a quench, and how quasimomentum is redistributed. By 
working in the weakly interacting limit we are able to 
produce analytic expressions. In particular, the time de- 
pendence of the static structure factor is quite simple. 
Our studies complement classical field calculations [12| 
and sophisticated numerically exact approaches, lar gely 
restricted to one-dimension [13l4l5| , or small systems [16| . 
Remarkably, much of the physics seen in the strongly in- 
teracting regime already appears in our calculations. For 
example, we find that the density correlations spread bal- 
listically [HI 12 1, and display damped oscillations. 

Formalism. — We consider a d-dimensional homoge- 
neous gas of bosons in an cubic optical lattice described 
by the single-band Bose- Hubbard Hamiltonian 18[ 3 : 
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where Oj(t) denotes the boson annihilation operator at 
site i, J denotes the hopping and U the on-site repul- 
sive interaction. The kinetic energy sum is over nearest 
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FIG. 1: (color online) Light-cone evolution of density- 
density correlations in ID. (c.f. Fig. 2, Ref. [l(|). Left: 
Density-density correlations function for a homogenous, non- 
interacting system <?];(t) = g d (t) (d — i — j) obtained from 
Eq. (0). Line shows the location of the peak in g d used 
to extract the velocity of spread of correlations. We find 
v = 3.7Ja, consistent with the spreading velocity expected 
for non-interacting particles. In the experiment of Ref. ffiol]) 
finds a velocity of v ~ 5.5Ja, in the strongly interacting limit. 
Right: Time-Evolution of the structure factor. Lighter colors 
indicate higher intensity. At t = 0, all momenta are equally 
occupied and S(k)(0) — 1 for all k. At intermediate times 
S(k) shows oscillations due to constructive interference of the 
atomic wave-packets when they travel a lattice unit. Higher 
momentum contributions to S(k) decay as 1/Jt, consistent 
with the linear spreading of correlations in real-space. 



neighbor pairs (ij). 

The basic objects of our study are the one and 
two body density matrices; <?j(t) = j(a] (t)a.j(t)) and 



al(t)a](t)ak(t)ai(t)), where ttj(t) denotes 



the boson annihilation operator at site i and time t. 
More generally we write the n body density matrix as 

atZ't^) = ^( a U t )-- a l( t ) a 3 J t )- a n( t ))- In various 
references, these are also referred to as the n body cor- 
relation functions, the 2n point functions, or the equal 
time Green's functions. 

The one and two body correlation functions arc read- 
ily measured in cold-atom experiments. The former is 
related to the momentum distribution function <?(k) = 

which i s p robed through 
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bandmapping |3|, H or time-of-flight The density- 

density correlation function g 1 ^ is related to the struc- 
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turc factor 5(q) = (p^P-q) = -Z) lJ eiq ' (l ^9% wnere 
Pq = 2k a k+q a k- The structure factor can be probed 
experimentally with Bragg scattering jH, 0] , or sophisti- 
cated imaging techniques [6|-|8[- 

The equations of motion for the n-body Green's func- 
tions are constructed from the equations of motion for 
the operators cii(t) and aj(t) : 



where all temporal dependence is implicit. 

For the one and two body Green's functions we obtain: 

id,g) = -J{g) +(J) - gf + Wisti - g%) (3) 
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where the notation (i) denotes a sum over all the near- 
est neighbors of site i in e?-dimcnsions. For example, 
in one dimension 9j + /j\ = 9j+i + 9j~i- ^ n a transla- 
tionally invariant system (such as the one we consider) 

9j+(j) ~ 9j ' anc ^ ^ e t erm proportional to J in Eq. (|3J) 
vanishes. 

The interaction term couples the n-body Green's func- 
tion with the n + 1-body Green's function. The full inter- 
acting many body dynamics is described by the resulting 
infinite set of coupled differential equations. 

Here we limit ourselves for the case of a shallow lat- 
tice, where interactions are weak following the quench. 
We first set (7 = and calculate the non-interacting 
density-density correlation functions (Eq. (0}). We then 
perturbatively include the effects of (7, determining how 
interactions influence quasi-momentum redistribution in 
the lattice (Eq. 0). 

Density- Density Correlations in ID. — We first con- 
sider a one-dimensional system as was studied in the ex- 
periment of Cheneau et al. [Io| . We choose a homoge- 
nous initial state with a density of hq bosons per site. 
At t = 0, the sites are completely decoupled, leading to 
a uniform quasi-momentum distribution with magnitude 
g{k) = n . 

In the absence of interactions, there is no quasi- 
momentum redistribution, and the momentum occupa- 
tions do not change in time. This can be easily seen by 
taking the Fourier transform Eq. ([3]). However density- 
density correlations given by Eq. (0]) show interesting dy- 
namics. 

To solve Eq. (|4]), it is convenient to define a reduced 
two body Green's function gV, — g^ — glgj — 9\ffL- As the 
two body correlations do not evolve in time, g 1 ^ obeys the 
same dynamical equation as gZ. At t = 0, g^i = n(n + l) 
if i = j = k = I, but otherwise. In Fourier space this 
becomes <?£f = uq{uq + l)6(p + q — r — s). 



Setting [7 = 0, Eq. (g]) is readily 
solved in Fourier space to yield g%j (t) = 

e -i2J t(cos(p)+cos(g)-cos(r)-cos(s))^pg^ t _ j n rea | 

space, the density-density correlation function g 1 ^ reads: 
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zeroth order Bessel function of first kind. 

In Fig. ([1}, we plot the dynamics of the two body 
Green's function as a function of time. As is apparent in 
the figure, the density correlations spread in a light-cone- 
like manner as observed in Ref. [ 1 01 ] . One can extract a 
characteristic velocity associated with the ballistic spread 
of correlations by plotting the location of the maximum 
of g d (indicated in Fig. Q] by the dashed line) as a func- 
tion of d. We obtain a velocity of v = 3.7 Ja. Note that 
there are other ways to extract a velocity from this data, 
which yield somewhat different values. 

Given that density correlations result from interfer- 
ence between left and right moving particles, one expects 
the correlations to spread at roughly twice the maxi- 
mum group velocity of a wavepacket v ~ 2i; max . Indeed, 



the dispersion e(k) = 



2v 



4Ja 



—2Jcos(ka) of this system yields 
The correlations spread about 10% 



slower than this, presumably due to contributions from 
other parts of the Brillioun zone. 

There is a striking similarity between the simple ana- 
lytic form derived in Eq. ([3]) for the density-density cor- 
relations in a non-interacting system and those observed 
experimentally and numerically in the strongly interact- 



ing regime 
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20J . This similarity is due to the fact that 
as long as the system is quenched sufficiently far from 
the critical point of the supcrfluid-insulator transition, it 
can be viewed as a dilute gas of weakly interacting quasi- 
particles, which propagate ballistically. One expects re- 
pulsive interactions to alter the propagation velocity of 
correlations, but have little effect on the functional form 
of the correlation function itself. 

In Fig. ([T]), we plot the structure factor obtained by 
taking the Fourier transform of the density-density cor- 
relation function. This can simply be read off from 
Eq. © as S(q)(t) = n (n + 1)J [4J tsin(g/2)] 2 . At 
t = 0, the structure factor is a constant as all momen- 
tum states are equally occupied. As the system begins 
to develop correlations between neighboring sites, the 
structure factor shows periodic oscillations whose ampli- 
tude decays in time. The oscillations indicate a rcphas- 
ing of the spreading wave-packets as they travel a lat- 
tice unit. Using the asymptotic behavior of the Bessel 
function Jq{z) ~ (2/nz) 1 / 2 cos(z — tt/4)) as z — >• oo, 
we find that for long times the oscillations have period 
T OS c = 7r/[4Jsin(fc/2)]. The envelope of of S(q) falls off 
as 1/Jt. 

A related experiment at MIT studied the Bragg in- 
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FIG. 2: (color online) Coherent redistribution of quasi- 
momentum to linear order in U/J in ID. Left: Redis- 
tribution of quasi-momentum (i%k = (atafc)) at times t = 
(thin, black), t = 0.25/ J (orange), 0.5/ J (red), 0.75 l/J 
(cyan), l/J (blue), 5/ J (purple) and 25/ J (solid, black) ob- 
tained by integrating Eq. ([6|. At short times, we find a coher- 
ent transfer of quasi-momentum from high momentum states 
to low momentum states. At longer times, we find a pile-up of 
particles near k = tt/2. Right: Spatial evolution of the one- 
body density matrix g](t) for different values of d = i — j. 
The correlations vanish if d is even. 



tensity of light reflected from a 3D cloud of expanding 
bosonic atoms released from an optical lattice [2j. Rel- 
atively small changes in their experiment would allow 
them to study the time evolution of the structure factor 
in a lattice. 

Momentum distribution in ID. — Following the quench 
we expect interactions to drive the system towards 
equilibrium. In particular, occupation of low momen- 
tum states should grow. To investigate this process 
we solve Eq. ([3]), replacing the two body correlator 
9ij(t) with the noninteracting result in Eq. J5]). Defin- 
ing x = 2Jt and using the integral identity J u (z) = 



^ f_ v d9e { ~ mf>+zcos ^ [2l|, we find that the occupation 
numbers obey: 



d t g q = Un (n + 1) i k J k {-x)J k {x) x (6) 

k— — OO 
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The right-hand-side of Eq. is invariant under the 
inversion q — > —q but switches sign under the transfor- 
mation q — > 7r — q. This implies that q = ±7r/2 is a sta- 
tionary point and states at q = n/2 have no dynamics. 
In real space, this symmetry implies = if d = \i — j\ 
is even. 

In Fig. @ we plot the evolution of the quasi- 
momentum states obtained by integrating Eq. ([6]). At 
t = 0, all momentum states are equally occupied. At 
short times following the quench, quasi-momentum states 
explore the band and the low momentum occupation 
begins to grow. Peaks develop near the stagnation 
points q = ±tt/2. The height of the peak scales as 



Uno(no + l)y/t/J. Due to this sharp spectral feature, in 
real space all (with d= odd) are appreciable at 

long times. 

The structure near q = ±7r/2 in Fig. [2] are reminis- 
cent of peaks seen in simulations of freely evolving ID 
hard core lattice bosons by Rigol and Muramatsu [lij . 
Despite the similarities, the relationship to our results 
is not clear, as we have no trap, and only weak inter- 
actions. Furthermore, they see a symmetric peak, while 
ours is antisymmetric. 

As shown in Fig. ©, at short times the single particle 
correlations spread in a manner similar to the density cor- 
relations. Local correlations are established on a time of 
order J -1 . Long range order, however, requires commu- 
nication between widely separated sites. Infinite range 
order is not found at any finite time. 

The linear approximation developed here breaks down 
whenever the height of the peaks near q = ±ir/2 be- 
come comparable to the initial density. At this point, 
we expect non-linear processes such as collisions be- 
tween quasi-momentum states at q and —q. This con- 
dition, (n q=7T /2(t) — no) /no ~ 1, occurs at time t ~ 
J/(U(no + l)) 2 . We do not know how robust the peaks 
are against collisions. One could argue that since two- 
body collisions do not change momentum states in ID, 
the quasi-condensate peaks may be long lived. 

Single site resolved experiments [l3| have studied the 
dynamics of the nearest neighbor correlation function fol- 
lowing a quench to shallow lattices (Vfinai < 2Er). They 
found that the nearest-neighbor single particles correla- 
tions scale linearly with U / J. 

Two dimensions. — We now generalize our results to 
higher dimensions. Concretely, we consider the case of 
a two-dimensional square lattice, initially containing no 
particles per site, and investigate the dynamics following 
a sudden reduction of the lattice depth to the weakly 
interacting limit. 

Repeating our previous ID arguments in two dimen- 
sions, we find that the structure factor is S(q x ,q y )(t) = 
n (n + l)Jo[4J t sin(q :E /2)] 2 Jo[4J t sm(q y /2)] 2 , which 
decays as 1/t 2 (as opposed to the 1/t decay in one- 
dimension) . This power law is a consequence of the bal- 
listic spread of correlations. After a time t, particle-hole 
correlations spread over a volume ~ u 2 t 2 where v is twice 
the characteristic velocity of an atom (v ~ 4y/2Ja). 

The time evolution of the momentum distribution plot- 
ted in Fig. ([3]) is quite different from the results in 
ID. Although the analog of Eq. © has four stationary 
points, {k x ,k y } = ±{7r/2, 7r/2}, these become zeros in 
Snk rather than peaks. At long times the distribution 
is characterized by a broad peak centered around k = 0. 
Thus to linear order in U , there is no condensation (char- 
acterized by a macroscopic occupation of the k = mode) 
nor quasi-condensation (characterized by a power law di- 
vergence of rife near k = 0). 

Conclusions and Outlook. — By considering the dy- 
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FIG. 3: Rapid equilibration of momentum distribution 
in two-dimensions Left: Momentum distribution along 
the {n, 7r} vector obtained by integrating Eq. © in two- 
dimensions assuming an initially uniform distribution. In 
contrast to the one-dimensional case, the distribution evolves 
rapidly to a broad peak at k = 0, with no further dynamics. 



namics of lattice bosons following a quench to a weakly 
interacting final state, we have explored how correlations 
develop in a many-body system. Remarkably we find 
that much of the behavior seen in the strongly interacting 
system is already present for weak interactions. For ex- 
ample, we find temporal oscillations and decay in the den- 
sity correlations, similar to those seen in experiments (loj . 
Within this framework we study how quasi-momentum 
states evolve following the quench, and investigate the 
initial stages of condensate formation. Our methods can 
readily be generalized to investigate fermionic systems, 
or more complicated lattice geometries, such as trian- 
gular and hexagonal lattices, that are at the frontier of 
cold-atom research [22[ ■ 

A major open question in the field of non-equilibrium 
dynamics is how off-diagonal long range order is gener- 
ated following a quench from a Mott insulating state. 
Over a decade ago there was a large body of work ask- 
ing analogous questions with thermal quenches 23j . The 
picture they developed was one of nucleation and sub- 
sequent coarsening. Similar physics is expected in the 
quantum case [24| . Experimental studies quenching into 
the regime h/U ~ hj J ~ 1ms find that in 3D a conden- 
sate develops on a timescale of order milliseconds [25[ . A 



key ingredient in the equilibration process is collisions, 
which first appear at order U 2 . 
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